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Abstract 

We present a normal form for travelling waves in one-dimensional excitable me- 
dia in form of a differential delay equation. The normal form is built around the 
well-known saddle-node bifurcation generically present in excitable media. Finite 
wavelength effects are captured by a delay. The normal form describes the be- 
haviour of single pulses in a periodic domain and also the richer behaviour of wave 
trains. The normal form exhibits a symmetry preserving Hopf bifurcation which 
may coalesce with the saddle-node in a Bogdanov-Takens point, and a symmetry 
breaking spatially inhomogeneous pitchfork bifurcation. We verify the existence of 
these bifurcations in numerical simulations. The parameters of the normal form are 
determined and its predictions are tested against numerical simulations of partial 
differential equation models of excitable media with good agreement. 
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Excitable media are often found in biological and chemical systems. Exam- 
ples of excitable media include electrical waves in cardiac and nerval tissue 
I jli. |2j|, cAMP waves in slime mold aggregation [3j and intracellular calcium 
waves [4j. Excitable media support localized pulses and periodic wave trains. 
In 2 dimensions rotating vortices (or spirals) and in 3 dimensions scroll waves 
[H ini IZl m ISj are possible. The critical behaviour of pulses, wave trains and 
spirals, i.e. propagation failure, is often associated with clinical situations. 
The study of spiral waves is particularly important as they are believed to 
be responsible for pathological cardiac arrhythmias [lOj. Spiral waves may 
be created in the heart through inhomogeneities in the cardiac tissue. Some 
aspects of spiral wave break up can be studied by looking at a one-dimensional 
slice of a spiral i.e. at a one-dimensional wave train [llj. 

We investigate critical behaviour relating to one-dimensional wave propaga- 
tion. We develop a normal form which allows us to study the bifurcation 
behaviour of critical waves. In particular, the normal form predicts a Hopf 
bifurcation and a symmetry breaking pitchfork bifurcation. The symmetry 
breaking pitchfork bifurcation can be numerically observed as an instability 
where every second pulse of a wave train dies. This seems to be related to 
alternans jl2|, ITB J , which are discussed in the context of cardiac electric pulse 
propagation. 

1 Introduction 

Many chemical and biological systems exhibit excitability. In small (zero-dimensional) 
geometry they show threshold behaviour, i.e. small perturbations immediately decay, 
whereas sufficiently large perturbations decay only after a large excursion. This behaviour 
is crucial for the electrical activation of cardiac tissue or the propagation of nerve pulses 
where activation should only be possible after a sufficiently large stimulus. Moreover, the 
decay to the rest state allows for the medium to be repeatedly activated - also crucial 
for the physiological functioning of the heart and the nervous system. One-dimensional 
excitable media support travelling pulses, or rather, periodic wave trains ranging in wave- 
length L from the localized limit L — »• cx) to a minimal value L^. below which propagation 
fails. Pulses and wave trains are best-known from nerve propagation along axons. In two 
dimensions one typically observes spiral waves. Spirals have been observed for example 
in the auto-catalytic Belousov-Zhabotinsky reaction [S], in the aggregation of the slime 
mold dictyostelium discoideum [3] and in cardiac tissue |2]. 

For certain system parameters the propagation of isolated pulses and wave trains may 
fail (see for example OHHI)- The analytical tools employed to describe these phenom- 
ena range from kinematic theory E|; asymptotic perturbation theory fSl ^1 120] 
to dynamical systems approaches [211 1221 123 • Numerical observations reveal that inde- 
pendent of the detailed structure of a particular excitable media model, the bifurcation 
behaviour of excitable media is generic. For example, a saddle-node bifurcation is generic 
for single pulses and for wave trains. However, there exists no general theory which ac- 
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counts for all bifurcations which may appear. In this paper we develop a normal form for 
excitable media which is build around the observation that the propagation failure of a 
one-dimensional wave train is mediated by the interaction of a pulse with the inhibitor of 
the preceding pulse. 

In Section 121 we briefly review some basic properties of excitable media and illustrate them 
with a specific example. In Section El we introduce the normal form. The properties and 
the bifurcation scenarios of this normal form are investigated in Section |31 In Section 
we show how the parameters of the normal form can be determined from numerical sim- 
ulations of the excitable medium and then compare the predictions of the normal form 
with actual numerical simulations of a partial differential equation. The paper concludes 
with a discussion in Section IHl 

2 One-dimensional excitable media 

Most theoretical investigations of excitable media are based on coupled reaction-diffusion 
models. We follow this tradition and investigate a two-component excitable medium with 
an activator u and a non-diffusive inhibitor v described by 

dtu = Duxx + ^{u,v), J-'{u,v) = u{l — u){u — Us — v) 

dtv = e {u — a v) . (1) 

This is a reparametrized version of a model introduced by Barkley [23 • Note that the 
diffusion constant D is not a relevant parameter as it can be scaled out by rescaling 
length. Although the normal form which we will introduce in Section |H1 is independent 
of the particular model used, we illustrate some basic properties of excitable media using 
the particular model (^. Later in Section Owe show correspondence of the predictions 
of our normal form with numerical simulations of the model (P). Our choice of model is 
motivated by the fact that this model incorporates the ingredients of an excitable system 
in a compact and lucid way. Thus, for > the rest state Mq = f o = is linearly stable 
with decay rates ai = Us along the activator direction and (J2 = ea along the inhibitor 
direction. Perturbing u above the threshold Us (in OD) will lead to growth of u. In the 
absence of v the activator would saturate at u = 1 leading to a bistable system. A positive 
inhibitor growth factor e and a > forces the activator to decay back to u = 0. Finally 
also the inhibitor with the refractory time constant (e a)~^ will decay back to f = 0. For 
a> 1/(1 — Us) the system is in zero-dimensional systems no longer excitable but instead 
bistable. 

In order to study pulse propagation in one- dimensional excitable media it is useful to 
first consider the case of constant v. The resulting bistable model is exactly solvable |2H1 
and the pulse velocity is Cf{y) = \f^[l — 2{us + v)]. Hence, excitability requires that Ug is 
below the stall value |. The quantity A = \ —Us characterizes the strength of excitability 
and c/(0) coincides with the solitary pulse velocity for e ^ 0. 

Clearly, for Ug < Uc = ^ and not too large a, pulse propagation fails for e larger 
than some ec- The critical growth factor Cc marks the onset of a saddle-node bifurcation 
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Figure 1: Typical profiles of the activator u (continuous line) and inhibitor v of close- 
to-critical travelling wave solutions in a one-dimensional ring of (^. Parameters are 
a = 0.22, Us = 0.1, e = 0.035. (a) Single pulse with L = 137.5. (b) Wave train with 
L = 550. 



|16|ITH1 I^. The saddle- node can be intuitively understood when we consider the activator 
pulse as a heat source, not unlike a fire-front in a bushfire. Due to the inhibitor the width 
of the pulse decreases with increasing e. Hence, the heat contained within the pulse 
decreases. At a critical width, or critical e, the heat contained within the pulse is too 
small to ignite/excite the medium in front of the pulse. 

Even if a given set of equation parameters allows for propagation of a single isolated 
pulse, the system may not necessarily support a pulse in a periodic box of finite length 
or a wave train consisting of several of such pulses: If the distance L between two con- 
secutive pulses of the train becomes too small, the pulses run into the refractory tail of 
the preceding pulse (see Figure E)), and may consecutively decay . Hence, propagation 
failure for periodic wave trains is controlled by the decay of the inhibitor, and propagation 
is only possible when the inter-pulse distance L is larger than a critical wavelength Lc. 
Note that Lc diverges for a —>■ when the decay rate of the inhibitor o"2 vanishes. The 
critical wavelength is a lower bound for the wavelength for the existence of periodic 
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wave trains. One can also think of keeping L fixed and, as before, vary e. Then the 
saddle-node e^L) is a monotonically increasing function. 

In the next Section we present a normal form which incorporates the saddle-node bifur- 
cation, and moreover predicts other types of bifurcations. 

3 The normal form 

It is well known that an isolated pulse undergoes a saddle-node bifurcation above a certain 
threshold value of the refractoriness ec- The corresponding generic normal form for such 
a saddle-node bifurcation is given by 

dtX = -fi-gX\ (2) 

where X is, for example, the amplitude or velocity of a pulse with the corresponding 
values at the saddle-node bifurcation subtracted. The bifurcation parameter /i measures 
the distance from the bifurcation point and is proportional to e — ec- The normal form Q 
accurately describes the behaviour of isolated solitary pulses in one- dimensional excitable 
media close to criticality. We use this normal form for a saddle-node bifurcation for an 
isolated pulse as a seed to construct a normal form for general travelling in excitable 
media incorporating finite wavelength effects. 

In particular we look at a single pulse on a ring with finite length, i.e. in a one- 
dimensional box with periodic boundary conditions, and at wave trains with a finite 
wavelength. Both cases are depicted in Figure |21 In that case the saddle-node will be 
disturbed and will depend on the length of the periodic box in the case of a single pulse 
(Figure Eti-) or on the wavelength of the wave train (Figure Eb). The interaction of the 
pulse (or a member of a wave train) with the preceding pulse (or more accurate with its 
inhibitor; see Figure EI) modifies as discussed above the bifurcation behaviour. We may 
therefore extend the saddle-node normal form to 

dtX = -fi- gX^ - f3oV{t - t) , 

where V{t — r) describes the inhibitor of the preceding pulse who is temporally displaced 
by r = L/cq where L is the wavelength of the pulse train, i.e. the distance between two 
consecutive pulses, and cq is its uniform velocity. We neglect here a possible temporal 
dependency of r. Note that in the case of a single pulse in a ring, V{t — r) describes the 
inhibitor of the single pulse which had been created by the pulse at the time of the last 
revolution around the ring and L is simply the length of the periodic box. 
We assume an exponential decay (in space and time) of the inhibitor of well separated 
pulses. This is the case for the system dT}. We may write Vit — r) = exp(— /cr)/i(X(t — r)) 
where the function h and the decay-rate k depend on the particular model chosen; for 
example for the model (0) we have k = ea. In the limiting case of isolated pulses we 
note that r ^ oo and V(t — t) 0, and hence we retrieve the unperturbed saddle- node 
bifurcation The ansatz for V{t — T) is a simplification where we ignore the cumulative 
effect of the inhibitor. (Note that the equation for the inhibitor f in ([Q) can be solved 
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directly and involves an integral over u, i.e. involves 'history'.) The unknown function 
h{X(t — r)) can be Taylor-expanded around the saddle-node X = 0. 
We summarize and arrive at the following normal form 

^^X = -^l-gX'-f3{J + X{t-T)), 

where (3 = /3oexp(— A;r) with k = ea for the model equation ([[]). Pulses have a nonzero 
width u which implies that the temporal delay r = L/cq has to be modified to r = 
{L — z/)/co. This equation already produces qualitatively all the results we will present 
in the subsequent Sections. However, much better quantitative agreement is achieved by 
taking into account that a variation in the amplitude implies a change in velocities and 
henceforth a change of the effective inhibition. If we allow for a single pulse to have a 
temporarily varying pulse amplitude or, in the case of a wave train consisting of distinct 
members, if we allow for different amplitudes of individual members of the wave train, we 
have to take into account that the propagation behaviour is amplitude dependent: Larger 
pulses have larger velocities. Hence, a pulse X{t) which is larger than its predecessor 
X{t — t) runs further into the inhibitor-populated space created by its predecessor. If 
X(t — r) < X{t) the finite wavelength induced shift of the bifurcation is stronger compared 
to the case of equal amplitudes. Conversely, if X(t — r) > X{t) the finite wavelength 
induced shift of the bifurcation is weaker compared to the case of equal amplitudes. 
This effect is stronger the larger the difference of the two amplitudes X{t — t) — X{t). 
The inclusion of the amplitude differences affects the bifurcation behaviour depending 
continuously on the difference X{t — t) — X{t). We thus add a term 7i(X(t) — X{t — t)) 
with < 7i <^ 1 into the wavelength dependant inhibitor term in (jS)), and arrive at 

dtX = -/i - gX^ - /3(7 + (1 - 7i)X(t - r) + ^^X{t)) , 

or after relabeling of /?, 7, 71 

dtX = -fi- gX^ - /5(7 + X{t - r) + -f,X{t)) . (3) 

It is this equation which we propose as a normal form to study bifurcations of one- 
dimensional wave trains. 

4 Properties of the normal form 

Before we show how to determine the parameters of the normal form, we will describe 
its properties with a main emphasis on bifurcations. Besides the well-known saddle-node 
bifurcation we identify a symmetry preserving Hopf bifurcation and a symmetry breaking 
spatially inhomogeneous pitchfork bifurcation. Numerical integration of partial differen- 
tial equation models of excitable media such as (0) confirm these bifurcation scenarios of 
the normal form Q. Although some of these bifurcations have been previously observed 
in numerical simulations, up to now there did not exist a unified framework to study 
these bifurcations. The normal form is able to identify these bifurcations as being generic 
for excitable media, rather than as being particular to certain models of excitable media. 
This is the main achievement of our present work. 



5 



4000 ^. SOOO 12000 



Figure 2: Temporal behaviour of the maximal amplitude U of the activator u for model ^ 
just above the subcritical Hopf bifurcation. The parameters are e = 0.3755 and L = 246. 
The other parameters are as in Figure |21 The subcritical character of the Hopf bifurcation 
is clearly seen at this length. If the length L is chosen to be slightly larger and closer to 
the actual bifurcation point, the oscillations appear to be stable for some time as shown 
in the inset for L = 246.9695 > Lc- However, the maximal amplitude U for L = 246.9695 
will eventually become visibly unstable and the pulse will die resulting in U = 0. 



4.1 Saddle- node bifurcation 

Numerical simulations of excitable media show that the bifurcations of a single propagat- 
ing pulse in a ring (as in Figure |2K) are different from the bifurcations of a wave train 
consisting of several distinct pulses (as in Figure |2b)- We first look at a single propagating 
pulse before in Section we look at the interaction of different pulses in a wave train. 
Equation has the following stationary solutions 

Xi,2 = ^M(l +7i) ± v'/3'(l + 7i)'-4^7(/i + /37)] . (4) 

The upper solution branch is stable whereas the lower one is unstable. The two solutions 
coalesce in a saddle-node bifurcation with 

^5iv = -^(l + 7i) at /i5iv = ^(l + 7i)'-/37- (5) 

Since j3 = /3oexp(— ear) is small we have fisN < 0. This indicates that the saddle-node 
of a periodic wave train occurs at smaller values of the bifurcation parameter fi than for 
the isolated pulse, and the bifurcation is shifted to the left with respect to the isolated 
pulse (see Figure . This is a well known fact which we numerically verified. Note that 
the limiting case of an isolated pulses with L oo implies r = and hence, /? = 0. The 
saddle-node of the isolated pulse with Xsn = at /i = described by (j2I) is recovered. 

Besides this stationary instability the normal form (jH)) also allows for a non-stationary 
bifurcation which we investigate in the next Section. 
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4.2 Symmetry preserving Hopf bifurcation 

The stability of the homogeneous solution X with respect to small perturbations of the 
form 6X expat can be studied by linearizing the normal form around X. We obtain 

a + 2gX + /57i + /3e~"" = . (6) 

Besides the stationary saddle-node bifurcation at a = (cf. (0)), a Hopf bifurcation 
(X = iuj is possible with 

u = PsinuT (7) 

Xhh = -^(coscur + 7i) . (8) 

In anticipation of the study of wave trains consisting of several distinct pulses we call this 
Hopf bifurcation of a single pulse in a ring a symmetry preserving Hopf bifurcation. From 
((Tj) we infer that a Hopf bifurcation is only possible provided /3r > 1, i.e. if the coupling 
is strong enough and the pulse feels the presence of the inhibitor of the preceding pulse 
sufficiently strong. Since Xhh > ^sn the symmetry preserving Hopf bifurcation sets in 
before the saddle-node bifurcation, independent of the value of p. Moreover, the Hopf 
bifurcation branches off the upper stable branch of the homogeneous stationary solutions 
(0}. In Figure we show a schematic bifurcation diagram with the saddle-node bifur- 
cation and the subcritical Hopf bifurcation for a single pulse in a ring. 

In numerical simulations of the Barkley model and also the Fitzhugh-Nagumo equations 
we could verify this scenario for a single pulse in a ring. A Hopf bifurcation had 
been previously observed numerically for the Barkley model |21] and in for the 
modified Barkley model (^. Hopf bifurcations have also been reported to occur in several 
other models of excitable media. In |2H1 1211 1101 a Hopf bifurcation was found in the 
8- variable Beeler- Renter model |3T] , and in [221 in the 4- variable Noble model [23] and 
in the 2- variable Karma model . We show here that Hopf bifurcations are generic for 
travelling waves in excitable media. 

In numerical simulations of model (Q) the Hopf bifurcation was found to be subcritical. 
Typical temporal behaviour of the maximal amplitude of the activator for model (0) close 
to the bifurcation is shown in Figure I^TI The inset shows the maximal amplitude slightly 
above the bifurcation point for about 20 periods. We counted more than 500 periods 
before stability was visibly lost and the maximal amplitude collapses to zero. We note 
that this may have easily lead to the wrong conclusion that the bifurcation is in fact 
supercritical rather than subcritical. 

4.3 Bogdanov-Takens point 

The saddle-node and the symmetry preserving Hopf bifurcation coalesce in a co-dimension- 
2 Bogdanov-Takens bifurcation for ujt 0. At the Bogdanov-Takens point we have 
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Figure 3: Sketch of the bifurcation diagram for a single pulse in a ring showing a stationary 
saddle-node bifurcation (SN) and a subcritical symmetry preserving Hopf bifurcation 
(HH). 



(3t = 1. The Hopf bifurcation and the saddle- node bifurcation have been suggested be- 
fore to be an unfolding of a Bogdanov-Takens point in |27j and later in [22] • The normal 
form provides a framework to study this unfolding. We were able to numerically verify 
the condition /5r = 1 derived from (|7j) by simulating the full partial differential equation 
(QJ. The parameters f3 and r will be determined further down in Sections El We have also 
numerically simulated the Fitzhugh-Nagumo equations to check that this bifurcation 
is not particular to our chosen model (^. 

The Bogdanov-Takens point is apparent in our normal form Q and can be derived from 
it. Close to the saddle- node and the Hopf bifurcation when ujt 0, the dynamics 
exhibits critical slowing down. We may therefore expand X{t — t) = X{t) — TdtX{t) + 
{t"^ /2)dttX{t) + 0{t^). The normal form (jS)) becomes at the Bogdanov-Takens point 

dtX = Y 

dtY = -aY -hX-^X\ (9) 

where X = X — Xi and Xi satisfies the stationary version of the normal form (jHl) and 
is given by The linear part of exhibits the correct eigenvalue structure of a 

Bogdanov-Takens point. The bifurcation parameters, a = 2(1 — /3r)/(r^/3) and b = 
2{P{l+'yi)+2gXi)/ (t'^P), measure the distance from the Hopf bifurcation and the distance 
from the saddle-node bifurcation, respectively. 

4.4 Spatially inhomogeneous pitchfork bifurcation 

Numerical simulations of systems such as (P) reveal that a group of several pulses in a ring 
do not undergo a sjTiimetry preserving Hopf bifurcation on increasing the refractoriness 
e, but instead develop a symmetry breaking, spatially inhomogeneous instability whereby 
every second pulse dies. In Figure 14.41 we show an example of such an inhomogeneous 
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Figure 4: Space-time plot of the activator u{x,t) demonstrating the spatially inhomoge- 
neous pitchfork bifurcation. Equation parameters are e = 0.04897125 and L = 300; all 
other parameters are as in Figure |21 Initially there is an initial group of 8 pulses which 
after some time exhibits a spatial period-doubling instability and subsequently evolves 
into a stable propagating wave train consisting of 4 remaining pulses. 



instability. Spatially inhomogeneous bifurcations have been observed before for periodi- 
cally paced excitable media [211 ESI EEl EZj • Here we show that this bifurcation is generic 
for wave trains in excitable media and does not require external pacing. 

The above mentioned inhomogeneous alternating instability is contained in our normal 
form. To investigate spatial instabilities we need to distinguish between consecutive 
pulses. We may rewrite the normal form as 

dtXi = -fi- gXf - /5(7 + Xi_,{t - r) + iiXi{t)) , (10) 

where the subscript / numbers the pulses in a wave train which interact with their nearest 
neighbours. Linearizing around the homogeneous solution Xi = X according to 

Xi=X + Se^'e'P^ and Xi^i = X + Se'^'e'^^^'^^ 

yields as a condition for stationary instabilities (ie. a = 0) 

X = — ^(71 + cosp — i sinp) . 

Hence stationary instabilities are possible for p = and for p = tt. In the homogeneous 
case p = 0, the instability is yet again the spatially homogeneous saddle-node (0). For 
p = 71 this is a new type of instability, and we have X = /5(1 — 7i)/2(yf at the bifurca- 
tion point. This spatially inhomogeneous bifurcation will be identified further down as a 
subcritical pitchfork bifurcation. The criterion p = vr for the inhomogeneous bifurcation 
is corroborated by numerical simulations where every second pulse dies within a wave 
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train (see Figure Note that in traveling wave coordinates of a partial differential 

equation model for excitable media, this instability would correspond to a subcritical 
period-doubling bifurcation. 



In order to study this p = vr-bifurcation within our normal form we need to consider two 
populations of pulses (/ = 1,2), X and which interact via their inhibitors with each 
other. We extend our normal form for the case p = tt to 

dtX = -^-gX^-(3{^ + Y{t-T)+^^X{t)) 

dtY = -^-gY^ -l3{^ + X{t-T)+^{Y{t)) . (11) 

The system flll|) for wave trains supports two types of stationary solutions; firstly the ho- 
mogeneous solution (jH), Xh = Yh, which may undergo a saddle-node bifurcation described 
by ®- There exists another stationary solution, an alternating mode, with 

Xa = -Ya + ^{l-^l) . (12) 

Associated with this solution is a pitchfork bifurcation at 

fJ'PF = 7 P7 = I^SN < fJ'SN , (13) 

4^ 9 9 9 

when 

Xpp = YpF=^{l-^i) . (14) 

Comparing ^ with ()1H|1 shows that the pitchfork bifurcation sets in before the saddle- 
node bifurcation. The upper branch of the homogeneous solution Xh given by (jH) at 
the pitchfork bifurcation point fipp coincides with fll4p . Hence the pitchfork bifurcation 
branches off the upper branch of the homogeneous solution. It is readily seen that the 
pitchfork bifurcation is subcritical because there are no solutions Xa possible for fi > fipp. 

We now look at the stability of the homogeneous solution X = Y = X^ = We 
study perturbations X = X^ + x exp at and Y = X^ + y exp at. Linearization yields as a 
condition for nontrivial solutions x and y. 

{a + 2gXt, + (3^,) = ±(3e-''^ . (15) 

The upper sign refers to an antisymmetric mode x = —y whereas the lower sign refers to 
a symmetric mode x = y. Stationary bifurcations occur at a = 0. The symmetric mode 
then coincides with the saddle-node bifurcation (0) whereas the antisymmetric mode ter- 
minates at the pitchfork bifurcation ()14j) . 

Non-stationary Hopf bifurcations are possible if a = iu. We then have 

uj = =F/5 sincjT (16) 



10 



X 




PF 



2g 



(1-Yl) 



PF 



sn; 



SN 



P 

2g 



(1+Yi) 



Figure 5: Sketch of the bifurcation diagram for a wave train in a ring showing a stationary 
saddle-node bifurcation (SN) and a subcritical pitchfork bifurcation (PF). 



One has physical solutions with a single-valued positive uj only for the symmetric case 
(the lower signs) which reproduces our results ((Zj) and (jH)) for the symmetry preserving 
Hopf bifurcation. For cur — the Hopf bifurcation moves towards the saddle-node (0) 
and coalesces with it at /3r = 1 in a Bogdanov-Takens point as described in Section |3] 
For lOT — > vr the limiting value of is = /?(1 — 7i)/(2(yf) which coincides with the 
pitchfork bifurcation Xpp in another codimension-2 bifurcation. At this bifurcation the 
Hopf bifurcation has a period T = 2r which corresponds exactly to the inhomogeneous 
pitchfork bifurcation with p = tt whereby every second pulse dies. 

For values ujt G [0, vr) the Hopf bifurcation always comes after the pitchfork bifurca- 
tion which has been numerically verified with simulations of the full system (^. 

This allows us to sketch the full bifurcation scenario for a wave train in a periodic ring as 
depicted in Figure lOl 

5 Determination of the parameters of the normal 
form 

In this Section we determine the parameters of the normal form from numerical sim- 
ulations of the full partial differential equations ([Q). We determine the free parameters 
z/, n, g, 7, 7i and f3o- We are then in the position to test how well the normal form Q 
reproduces the solution behaviour of the full partial differential equation (jT)). We use here 
as equation parameters for ([T]) a = 0.22, = 0.1 and D = 1. 

The parameter u which modifies the delay time r due to the finite width of a pulse is 
easily determined as a typical width of a pulse in the parameter region of interest. We 



and 



Xh = —(± cos u;r - 71) . 



(17) 
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find z/ = 32. We note that there is some ambiguity in the determination of z/, and one 
may as easily justify v G [29, 34]. 



The parameters \i and g can be determined by studying the isolated pulse with L oo. 
The normal form reduces to 

dtX = gX^ . (18) 

We have [x = a{e — ec) with ec being the critical e at the saddle-node. Solutions of (fTHj) 
are obtained by quadrature 



X{t) 




tanh{^y—fig{t — to)) for < 



X{t) = -J-tan{^{t - to)) for /x > . 



(19) 



For small deviations from the saddle-node X = this solution may be expanded to obtain 
X{t) /it which obviously corresponds to the solution of equation ()18|) linearized around 
the saddle-node. The solution (fT9|l has its inflection point at the saddle-node X = 
where its slope is /i. We can therefore determine fi by measuring dtX at the inflection 
point for different values of (e — ec). This allows us to determine a via a = fi / {e — ec) . In 
Figure El we show the results of dtX versus (e — ec). The numerical results are obtained 
by letting a stable pulse which was created at some e < ec, decay in an environment with 
e > ec. The relaxation then allows us to determine the slope at the saddle-node. Using a 
least-square fit we obtain a = 1.455. 



The parameter g can now be determined by looking at the stationary problem dfX = 0. 
The behaviour of the amplitude of the activator close to the saddle-node versus e is de- 
picted in Figure El It clearly demonstrates quadratic behaviour typical for saddle-nodes. 
The normal form for the saddle-node of an isolated pulse Q yields (e — ec) = {g/a)X'^ 
which we can use upon using the above measured value of a to obtain g = 0.31 from a 
least square fit. 

To determine the missing parameters /3o, 7 and 71 we need to study a pulse in a periodic 
ring of finite length L. The self-interaction of the pulse with its own inhibitor modifies 
the saddle-node bifurcation as discussed in Section |21 We will look here at the symmetry 
preserving Hopf bifurcation. To study the Hopf bifurcation we study the circulation of a 
single pulse in a periodic ring instead of a wave train consisting of more than one pulse. In 
the latter case the subcritical pitchfork bifurcation sets in before the symmetry preserving 
Hopf bifurcation. 

Combining the expressions for the angular frequency u and the deviation Xhh of 
the pulse amplitude from the saddle- node at the bifurcation point, (|7|) and (jH]), we can 
eliminate the so far undetermined parameter jS to determine 71. We obtain 

— ^ = - — (cot(cjr) + 7i sin~^(cjr)) . (20) 
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Figure 6: The temporal derivative dtX at the saddle-node versus e — ec- The slope deter- 
mines the parameter a = /i/(e — Cc). The stars are obtained by numerically integrating 
the partial differential equation (^Q). The line is a least-square fit. 
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Figure 7: Bifurcation parameter e — ec versus amplitude X close to the saddle-node of the 
isolated pulse. The crosses are obtained by numerically integrating the partial differential 
equation (Q). The line is a quadratic fit. 



13 



"I 



4 - 




-2 



0.5 



1 . 5 



2 . 5 



LOT 



Figure 8: Plot of Xhh/^ versus uot. The crosses depict results from numerically integrat- 
ing the full partial differential equation (^. The line represents the analytical expression 
()2()j) of our normal form. This fit is achieved by choosing 71 = 0.31. 



In Figure El we show a plot of numerically obtained values for the quantity Xhh /^^ by 
integrating the partial differential equation and the result of our normal form (jSOl)- 
In the numerical simulations of (0) we measured the frequency of the subcritical Hopf 
bifurcation uj and X^h which is the difference between the amplitude at the symmetry 
preserving Hopf bifurcation and the saddle-node value of the isolated pulse. Expression 
(|2(Jj) matches the numerical simulations well for 71 = 0.31. 

We can now determine the parameters (3q and 7 by looking at the shift of the bifurca- 
tion parameter /i and the shift of the critical amplitude X at the saddle-node for finite 
wavelength L with respect to the values at the saddle-node for an isolated pulse with 
infinite wavelength. The saddle-node is shifted with respect to the isolated pulse with 
fiSN = XsN = 0. We express the shifted finite-L saddle-node by 



where 5 = S{L) expresses the shift in the bifurcation parameter at the saddle-node and 
^ = C,{L) expresses the shifted value of the amplitude when compared to the isolated 
pulse. 

The values for 6{L) and ^(L) can be measured by numerically solving ([Q) and determing 
the finite wavelength induced saddle-node. We did so by expressing (P) in travelling wave 
coordinates and treating the problem as a boundary value problem. The thereby obtained 
values for S{L) and C,{L) have to be compared with the corresponding expressions of the 
normal form. In the normal form (jH)) the finite wavelength induced shift is represented by 
the finite length correction 



^ + S = -giX-0\ 



r(t-r) =/?oexp(-ae(L-z/)/c)(7 + X(t-r)+7iX(t)) . 
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Figure 9: Shift of the bifurcation parameter (e — ec) at the saddle-node as a function of L. 
Crosses depict results from a numerical simulation of the full partial differential equation 
(0). The line represent the analytical expression of our normal form. Agreement of 
the two curves is achieved for 7/^0 = 0.115. 



Neglecting which is justified for not too large deviations from the isolated pulse, 

we obtain 



where v and 71 had already been determined. Note that the 0(^^)-term, g^'^ = /3^(1 + 
7i)/4(7 appears, of course, correctly in the shift of the bifurcation parameter f^sN in ®- 
Although the inclusion of 71 7^ is not necessary for the existence of a Hopf bifurcation 
(see ()7I8|) ). it is significant to obtain good quantitative agreement. Whereas for the sta- 
tionary bifurcations the inclusion of 71 is in effect a redefinition of (3, it is vital in the case 
of the Hopf bifurcation because it allows for a decoupled dependence of the frequency u 
and the amplitude Xhh on L. In Figure El we show a plot of the bifurcation parameter 
shift as a function of length L. Agreement of numerically obtained values from a sim- 
ulation of the full partial differential equation (0) with the expression derived from our 
normal form (j^H), which implies 6 = a{e — ec) is assured provided 7/^0 = 0.115. We recall 
that ec is the critical refractoriness at the saddle-node of an isolated pulse. In Figure we 
show a plot of the pulse amplitude shift as a function of length L. Agreement of numer- 
ically obtained values from a simulation of the full partial differential equation (P) with 
the expression derived from our normal form (j^^ is given provided (/?o/25')(l + 71) = 1.1. 
Combining these results we can solve for Pq and 7 and obtain 7 = 0.32 and Po = 0.53. 
This finalizes the determination of the free parameters of the normal form (jHl). 



S 



7/3o exp(— ea(L — z/)/c) 

^(l + 7i) exp(-ea(L - z/)/c) 



(21) 
(22) 
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Figure 10: Shift of the amphtude deviation X at the saddle-node as a function of L. 
Crosses depict results from a numerical simulation of the full partial differential equation 
(0). The line represent the analytical expression ()22|) of our normal form. Agreement of 
the two curves is achieved for {(3Q/2g){l + 71) = 0.91. 



We are now in the position to check the validity of our normal form by testing whether the 
normal form (jH} with the above determined parameters is able to reproduce observations 
of the numerical simulation of the full partial differential equation (^. We already noted 
in Section 0] qualitative agreement i.e. the correct bifurcation behaviour. Here we show 
quantitative agreement with the behaviour of (^. 

In particular we look at the symmetry preserving Hopf bifurcation described in Section HJ 
In Figure El and Figure El we show a comparison of the analytical results for the frequency 
and the amplitude at the Hopf bifurcation, ((7|) and (|HI), with results obtained from nu- 
merically integrating (0). The Figures show good agreement. Note that the parameters 
/So and v were not determined by fitting data representing the symmetry preserving Hopf 
bifurcation and henceforth the two figures, Figure Eland Figure El are indeed predictions. 

Figure El shows a comparison of numerical simulations of and our analytical results 
p4j) for the spatially inhomogeneous pitchfork bifurcation. For the determination of the 
parameters of the normal form we have not used any fitting which involved results from 
the spatially inhomogeneous pitchfork bifurcation. The agreement in Figure El therefore 
demonstrates that the normal form can indeed be used to obtain quantitative agreement 
and make quantitative predictions. 
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Figure 11: Amplitude deviation at the Hopf bifurcation X^h as a function of L. Crosses 
depict results from a numerical simulation of the full partial differential equation (^. The 
line represent the analytical expression (jS)) of our normal form. 
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Figure 12: Angular frequency at the Hopf bifurcation as a function of L. Crosses depict 
results from a numerical simulation of the full partial differential equation . Stars were 
obtained by solving the analytical expression ((7j) of our normal form for uj. 
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Figure 13: Amplitude deviations at the pitchfork bifurcation Xpp as a function of e. 
Crosses depict numerical results from integrating the full partial differential equation ((T}. 
Stars show the corresponding values calculated by using the normal form result p4|l . 



6 Summary and discussion 

We have constructed a normal form for travelling waves in one-dimensional excitable me- 
dia which takes the form of a delay differential equation. The construction is based on the 
well-known observation that the interaction of a pulse with the inhibitor of the preceding 
pulse modifies the generic saddle-node bifurcation of an isolated pulse. The normal form 
(jHl) exhibits a rich bifurcation behaviour which we could verify by numerically simulating 
the partial differential equation (0). Besides the well known saddle-node bifurcations for 
isolated pulses and for periodic wave trains the normal form also exhibits a symmetry 
preserving Hopf bifurcation and a symmetry breaking, spatially inhomogeneous pitchfork 
bifurcation. Moreover, the normal form shows that the saddle-node and the Hopf bifur- 
cation are an unfolding of a Bogdanov-Takens point as previously suggested in |2Z1I2S1- 
The symmetry preserving Hopf bifurcation is found to occur before the saddle-node bifur- 
cation for a single pulse in a ring. For a wave train consisting of several pulses in a ring, 
the Hopf- and the saddle-node bifurcations occur after the symmetry breaking pitchfork 
bifurcation in which every second pulse dies. We could verify these scenarios in numer- 
ical simulations of the modified Barkley-model (0) and the Fitzhugh-Nagumo equations 
PB] . These bifurcations have been observed before but had previously not been described 
within a unified framework of one normal form. 

We were able to determine the parameters of the normal form from numerical simulations 
of the partial differential equation (^. Using these numerically determined parameters 
we showed excellent agreement between the normal form and the full partial differential 
equation (^. We could quantitatively describe the symmetry preserving Hopf bifurcation 
and the inhomogeneous pitchfork bifurcation. Moreover, we were able to quantify the 
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Bogdanov-Takens bifurcation. 



The symmetry preserving Hopf bifurcation has been studied intensively before. It was 
observed numerically for example for the Barkley- model [27j. Interest has risen in the 
Hopf bifurcation in the context of cardiac dynamics because it leads to propagation fail- 
ure of a single pulse on a ring. It is believed to be related to a phenomenon in cardiac 
excitable media which goes under the name of alternans. Alternans describe the sce- 
nario whereby action potential durations are alternating periodically between short and 
long periods. The interest in alternans has risen as they are believed to trigger spiral 
wave breakup in cardiac tissue and ventricular fibrillation [T21 HH E21 UHl IHH] • Besides 
numerical investigations of the Barkley model |2Z] , the modified Barkley model , the 
Beeler-Reuter model [21111111211 EH], the Noble-model [2211201 and the Karma-model [S2!, 
where a Hopf bifurcation has been reported, there have been many theoretical attempts 
to quantify this bifurcation for a single-pulse on a ring. Since the pioneering work [T^ 
alternans have been related to a period-doubling bifurcation. It was proposed that the 
bifurcation can be described by a one- dimensional return map relating the action poten- 
tial duration [APD) to the previous recovery time, or diastolic interval (DI), which is 
the time between the end of a pulse to the next excitation. A period-doubling bifurca- 
tion was found if the slope of the so called restitution curve which relates the APD to 
the DI, exceeds one. A critical account on the predictive nature of the restitution curve 
for period-doubling bifurcations is given in [SHI IMl • In [21] the instability was analyzed 
by reducing the partial differential equation describing the excitable media to a discrete 
map via a reduction to a free-boundary problem. In the Hopf bifurcation could be 
described by means of a reduced set of ordinary-differential equations using a collective 
coordinate approach. In [TTl jSU] the bifurcation was linked to an instability of a single 
integro-delay equation. The condition for instability given by this approach states - as 
in some previous studies involving one-dimensional return maps - that the slope of the 
restitution curve needs to be greater than one. However, as evidenced in experiments [IHl 
and in theoretical studies (SHI EE] alternans do not necessarily occur when the slope of 
the restitution curve is greater than one. In further studies it would be interesting to see 
how our criterion /3r > 1 is related to that condition. Note that u is related to the APD 
and T = (L — z/)/c to the recovery time DI. 

In the context of alternans the Hopf bifurcation has been described as a supercritical 
bifurcation [1111321123120] (although their occurrence is related to wave break up [29j). 
Further study will explore whether the subcritical character of the Hopf bifurcation we 
find is model-dependant or indeed generic. 

To go beyond the case of a single pulse circulating in a ring, periodically stimulated 
excitable media have been studied in the context of alternans [HI |321 jSHl ISll jSHl 1^ - 
In [21 EOl one-dimensional maps were developed to study the Hopf bifurcation and the 
transition to conduction blocks. In [32] a nonlocal partial differential equation has been 
proposed to study spatiotemporal dynamics of alternans. It would be interesting for fur- 
ther studies to see how the transition to conduction blocks explored in these paced cardiac 
excitable media [311 ESI EHl E^ can be described by the spatially inhomogeneous pitch- 
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fork bifurcation we found in our normal form. The pitchfork bifurcation however does not 
require a fixed pacing site and does not require external pacing but rather is dynamically 
induced. This may aid in investigating the formation of conduction blocks purely as a 
dynamical phenomenon of wave trains. 

An interesting scenario in our normal is the coalescence of the spatially inhomogeneous 
pitchfork bifurcation with the Hopf bifurcation when ujt = vr. This condition implies 
T = 2r. Then the Hopf frequency is in resonance with the spatial instability in which 
every second pulse dies. Connections to alternans of this scenario are planned for further 
research. 

Ideally, one would like to deduce the normal form directly from a model for excitable me- 
dia and determine its parameters without relying on numerical simulations of a particular 
excitable medium. One initial path along that avenue could be to use the non-perturbative 
approach developed in [221 to determine the parameters. This method was developed to 
study critical wave propagation of single pulses and pulse trains in excitable media in 
one and two dimensions. It was based on the observation that close to the bifurcation 
point the pulse shape is approximately a bell-shaped function. Numerical simulations 
show that this is the case for the Barkley model ([Q) close to the saddle-node bifurcation. 
A test function approximation that optimises the two free parameters of a bell-shaped 
function, i.e. its amplitude and its width, allows us to find the actual bifurcation point, 
ec, and determine the pulse shape for close-to-critical pulses at excitabilities near ec- This 
method has so far also been successfully applied to other non-excitable reaction-diffusion 
equations jlHj |33] . To apply the method for our purpose is planned for future work. 
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